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INTRODUCTION 

This  study  investigated  a  novel  imaging  technology  of  digital  breast  tomosynthesis  (DBT)  to 
improve  early  breast  cancer  detection. 

Breast  cancer  is  a  major  problem  and  the  most  common  cancer  among  women.  Approximately 
over  10%  of  women  in  the  US  will  develop  breast  cancer  during  their  lifetime  L  2  and  30%  to 
40%  of  American  women  who  get  breast  cancer  die  from  it.  1  Annually,  a  total  of  348,000  cases 
of  breast  cancer  is  diagnosed  and  almost  115,000  are  killed  by  breast  cancer  in  the  US  and 
European  Community.  2  Mammography  is  currently  the  most  important  and  efficacious  tool  for 
the  early  detection  of  breast  cancer. 3  However,  the  nature  of  the  two-dimensional  mammography 
makes  it  very  difficult  to  distinguish  a  cancer  from  overlying  breast  tissues,  especially  for  dense 
breast  cases. 

Digital  breast  tomosynthesis  is  a  three-dimensional  breast  imaging  method  to  reconstruct  a  set  of 
planes  in  the  breast  from  limited-angle  series  of  projection  images  as  the  x-ray  source.  A  variety 
of  tomosynthesis  reconstruction  algorithms  have  been  proposed  including  the  traditional  shift- 
and-add  (SAA),  the  image-stretching  method  proposed  by  Nilklason  and  colleagues  4,  the 
maximum  likelihood  iterative  algorithm  (MLEM)  by  Wu  et  al.  5’  6,  tuned-aperture  computed 
tomography  (TACT)  reconstruction  methods  developed  by  Webber  and  investigated  by 
Suryanarayanan  et  al 7’ 8,  algebraic  reconstruction  techniques  (ART)  9~n ,  filtered  back  projection 
(FBP)  12~16,  and  matrix  inversion  tomosynthesis  (MITS)  16~18. 

The  purpose  of  this  project  is  to  optimize  and  compare  several  different  tomosynthesis 
reconstruction  algorithms  that  are  either  initially  investigated  in  our  lab  or  currently  very 
popular,  and  to  optimize  the  imaging  acquisition  techniques.  Based  on  this  project,  we  hope  to 
contribute  to  the  optimal  breast  tomosynthesis  technique  for  better  breast  cancer  detection. 

Note  that  although  this  predoctoral  fellowship  was  awarded  for  two  years,  it  has  now  been 
concluded  in  one  year  and  five  months.  The  recipient,  Ying  Chen,  just  graduated  in  August  2007 
and  has  resigned  from  the  fellowship.  The  fellowship  provided  her  with  a  solid  foundation  in 
cancer  research  and  has  allowed  her  to  continue  with  cancer  research  as  a  tenure-track  assistant 
professor  at  Southern  Illinois  University. 

BODY 

Task  1.  Optimization  of  different  candidate  tomosynthesis  reconstruction  methods  (Month 
1-12): 

1.1.  Select  candidate  algorithms  from  different  algorithm  families  for  optimization  and 
comparison.  2-4  candidate  algorithms  from  different  algorithm  families  will  be  chosen. 
(Month  1-4). 

This  task  has  been  completed  and  has  resulted  in  publications  (see  #1,  #2  and  #3  in  Reportable 
Outcomes).  We  have  investigated  a  few  reconstruction  algorithms,  including  shift-and-add 
(SAA)  algorithm,  Nilklason’s  image  stretching  shift-and-add  (NIKL),  maximum  likelihood 
expectation  maximization  (MLEM),  matrix  inversion  tomosynthesis  (MITS),  and  filtered  back 
projection  (FBP). 
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Among  those  algorithms,  MITS  was  originally  invented  in  our  lab  by  Dobbins  and  we  applied 
it  to  the  breast  tomosynthesis  imaging  successfully  !6’  17 .  We  also  developed  our  FBP 
reconstruction  based  on  central  slice  theorem  and  Fourier  frequency  sampling  density.  In  order  to 
control  the  high  frequency  noise  amplification,  Hamming  Gaussian  filters  were  designed  and 
applied  to  our  FBP  algorithm  16’  19 .  MLEM  algorithm  is  an  effective  iterative  method  in  breast 
tomosynthesis  reconstruction.  However,  it  is  time-consuming  due  to  intensive  computation. 
Therefore,  in  this  project,  we  selected  back  projection  (BP),  MITS  and  FBP  as  our  three 
candidate  algorithms  for  comparison  and  optimization.  But  we  also  did  related  research  on  other 
algorithms  such  as  SAA,  NIKL,  and  MLEM. 

During  our  investigation  of  different  algorithms,  we  found  that  quite  a  few  other  algorithms 
depend  on  a  traditional  SAA  method.  However,  a  simple  SAA  reconstruction  algorithm  is  not 
entirely  suitable  for  isocentric  motion  in  breast  tomosynthesis,  especially  for  small  structures 
such  as  microcalcifications,  which  have  an  important  bearing  in  clinical  breast  cancer  detection 
tasks.  A  point-by-point  back  projection  (BP)  is  important  for  isocentric  motion  in  breast 
tomosynthesis. 

We  have  published  our  results  in  Medical  Physics,  the  premiere  peer-reviewed  journal  in  the 
field  of  Medical  Physics;  please  see  Appendix  #1  for  the  reprinted  publication. 

1.2.  Characterize  the  effect  of  three  acquisition  parameters  including  total  Tomographic- 
Angle  (TA),  Number  of  projection  images  (N),  and  Reconstruction-Slice-Spacing  (RSS)  for 
each  reconstruction  algorithm,  according  to  physical  measurements  of  impulse  response 
analysis,  modulation  transfer  function  (MTF)  and  noise  power  spectrum  (NPS).  (Months  5- 
12). 

This  task  has  been  completed  and  has  resulted  in  publications  (see  #1,  #2,  #3,  and  #4)  in 
Reportable  Outcomes).  Parameters  of  a  selenium-based  direct  conversion  Siemens  Mammomat 
Novation  DR  prototype  system  was  modified  to  be  used  for  geometries  of  the  simulation.  Seven 
different  combinations  of  acquisition  parameters  were  investigated.  We  have  characterized 
candidate  algorithms  and  different  combinations  of  acquisition  parameters  by  simulation  and 
physical  measurements  of  impulse  response,  MTF  and  NPS. 

For  each  combination  of  acquisition  parameters,  datasets  of  tomosynthesis  projection  sequences 
with  several  different  impulse  locations  were  simulated  with  a  ray-tracing  method  and  then 
reconstructed  with  candidate  tomosynthesis  algorithms.  Impulse  response  was  analyzed  and 
compared.  We  have  measured  MTF  values  of  each  algorithm  and  acquisition  parameters  based 
on  two  parts:  1)  measured  system  MTF  of  the  detector  by  previously  published  edge -method  ' 
21 ;  2)  the  relative  reconstruction  MTF  associated  with  specific  reconstruction  algorithm  and 
acquisition  parameters.  For  MITS  and  FBP  algorithms,  due  to  applied  deblurring  algorithms, 
there  is  no  real  meaningful  DC  components  to  normalize  the  MTF.  Therefore,  we  calcualted 
relative  reconstruction  MTF  as  the  magnitude  of  Fourier  Transform  of  the  impulse  response 
associated  with  specific  reconstruction  algorithm  and  acquisition  parameters. 
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In  order  to  compute  the  noise  power  spectrum  (NPS)  for  specific  acquisition  parameters,  two 
identical  phantom  slabs  (47%  water /  53%  adipose  equivalent)  for  a  total  of  4  cm  thickness  were 
put  directly  on  the  surface  plate  (detector  cover)  to  mimic  the  breast  tissue  equivalent 
attenuation  and  scattered  radiation.  Ten  identical  tomosynthesis  sequences  of  flat  images  with  the 
phantom  slabs  on  the  detector  were  acquired.  The  tomosynthesis  sequences  were  then 
reconstructed  by  each  candidate  algorithm.  Mean-subtracted  NPS  was  analyzed  is  to  remove 

22  23 

fixed  pattern  noise,  including  structured  noise  and  system  artifacts.  ’ 

Figures  1(a),  1(b)  and  1(c)  show  an  example  of  in-plane  impulse  response  analysis  for  BP,  MITS 
and  FBP  algorithms  respectively  with  25  projection  images  and  ±25°  total  tomographic  angle  of 
x-ray  tube’s  movement.  The  impulse  was  located  25  mm  above  the  detector  cover  surface  plate 
and  near  the  chest  wall. 


0) 


E 


Figure  1(a):  Impulse  response  example  of  BP 


Figure  1(b):  Impulse  response  example  of  MITS 


Figure  1(c):  Impulse  response  example  of  FBP 

We  found  that  all  three  algorithms  can  reconstruct  the  simulated  impulses  to  give  sharp  in-plane 
responses.  MITS  and  FBP  have  tails  due  to  their  edge  enhancement  effects.  For  out-of-plane 
response,  compared  with  BP,  MITS  and  FBP  were  always  better  for  removal  of  out-of-plane  blur 
artifacts.  There  was  no  substantial  difference  for  BP  with  different  acquisition  parameters.  For 
MITS  and  FBP,  bigger  N  cases  were  better  in  high  frequency  noise  removal.  FBP  performed 
better  with  wider  tomographic  angle  of  tube’s  movement.  For  MITS,  when  impulse  is  located 
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near  chest  wall  and  with  a  narrower  tomographic  angle,  MITS  performed  well  with  less  artifacts 
associated  with  isocentric  motion  of  breast  tomosynthesis. 

Figures  2(a),  2(b)  and  2(c)  show  an  example  of  our  relative  MTF  analysis  of  BP,  MITS  and  FBP, 
respectively,  with  seven  different  combinations  of  acquisition  parameters.  Acquisition 
parameters  were  named  based  on  tomographic  angle  (TA)  and  number  of  projection  images  (N) 
for  clarity.  Scrub  (“Bo”)  means  a  simple  dark  current  image  between  each  frame.  Binning  mode 
(“Bin”)  means  to  sum  up  the  neighboring  2  pixel  values  along  the  x-ray  tube’s  motion  to  speed 
up  the  tomosynthesis  scan.  For  examples,  B0-bin-A50-P25  means:  Angular  range  =50°, 
Projection  images=25,  with  “Bo”  scrub  and  Binning  mode;  A44-P13  means:  Angular  range 
=44°,  Projection  images=13. 


BP:  all  techniques,  nearcw 


Figure  2(a):  Relative 
MTF  of  BP  for  several 
different  acquisition 
combinations 


MITS  Slabbing^:  all  techniques,  nearcw 
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Figure  2(b):  Relative 
MTF  of  MITS  for 
several  different 
acquisition 
combinations 
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Mean-subtracted  NFS  Mean-subtracted  NFS 


CorFBP  Hamm:  all  techniques,  nearcw 


Figure  2(c):  Relative 
MTF  of  FBP  for 
several  different 
acquisition 
combinations 


Figures  3(a),  3(b)  and  3(c)  show  an  example  of  our  mean-subtracted  analysis  of  BP,  MITS  and 
FBP,  respectively,  with  seven  different  combinations  of  acquisition  parameters. 


Figure  3(a):  Mean- 
subtracted  NPS  of  BP 
for  several  different 
acquisition 
combinations 


Figure  3(b):  Mean- 
subtracted  NPS  of 
MITS  for  several 
different  acquisition 
combinations 
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CorFBP  Hamm 


co 


Figure  3(c):  Mean- 
subtracted  NPS  of 
FBP  for  several 
different  acquisition 
combinations 


We  found  that  compared  with  other  acquisition  modes,  the  B0-bin-A50-P25  (N=25,  TA=±25° 
angular  range,  binning  mode)  has  lower  high  frequency  noise  characteristics  for  BP  and  MITS. 
Among  the  other  six  non-binning  acquisition  modes,  BP  performed  slightly  worse  for  N=13 
cases.  MITS  performed  worse  for  less  number  of  projection  images  and  wider  tomographic  angle 
of  N=13,  TA=±22.5°  case;  MITS  had  better  noise  characteristics  for  bigger  number  of  projection 
images  and  wider  tomographic  angle  of  N=49,  TA=±25°  case.  With  same  acquisition  parameters, 
the  SAA  and  FBP  showed  no  difference  when  we  varied  the  RSS.  For  MITS  reconstruction, 
compared  with  2mm  and  4mm  slice  spacing,  1mm  slice  spacing  has  better  noise  response  at  low- 
to-middle  frequency  range.  At  middle  to  high  frequency  range,  their  noise  responses  were 
similar. 

We  have  proposed  an  image  reconstruction  algorithm  of  Gaussian  Frequency  Blending  (GFB) 
algorithm  using  NPS  analysis  to  combine  MITS  and  FBT  together  for  better  image 
reconstruction  of  breast  tomosynthesis.  We  are  working  on  a  Medical  Physics  journal  manuscript 
preparation  for  GFB  algorithm. 

We  have  used  impulse  response  and  MTF  analysis  method  to  compare  BP  and  traditional  SAA 
algorithms  (see  Appendix  #1,  #2,  and  #3).  Our  results  of  impulse  response  analysis  for  BP  with 
relevance  to  morphology  of  structures  such  as  microcalcifications  were  published  in  Medical 
Physics,  the  premiere  peer-reviewed  journal  in  the  field  of  Medical  Physics,  please  see  Appendix 
#1  for  the  reprinted  publication. 

Task  2.  Comparison  of  different  candidate  tomosynthesis  reconstruction  methods  (Months 
13-24): 

To  allow  for  computational  modeling,  human  observer  study  efforts  were  scheduled  for  the 
second  budget  year.  However,  we  are  relinquishing  seven  months  of  the  second  year  because  the 
predoctoral  fellowship  recipient,  Ying  Chen,  has  just  received  her  PhD  and  graduated. 

Instead,  we  have  published  a  methodology  of  relative  NEQ(f)  analysis  to  compare  the  efficacy  of 
different  candidate  tomosynthesis  reconstruction  methods  and  have  completed  the  comparison 
tasks.  The  methodology  of  NEQ(f)  analysis  was  published  as  a  proceedings  paper  at  SPIE  (see 
Appendix  #4).  It  combines  MTF  of  signal  performance  and  NPS  of  noise  characteristics  to 
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enable  one  to  evaluate  the  performance  of  different  acquisition  parameters  and  algorithms  for 
comparison  and  optimization  purposes.  For  each  acquisition  parameter  and  candidate 
reconstruction  algorithm,  the  same  inputs  were  used  in  all  cases  (the  same  simulated  impulse 
magnitude  for  relative  reconstruction  MTF  measurement  and  same  accumulative  tomosynthesis 
sequence  exposure  level  for  NPS  measurement).  Therefore,  one  can  make  relative  comparisons 
to  evaluate  performance  of  different  algorithms  and  acquisition  modes  with  the  same  inputs. 
Figures  4(a),  4(b)  and  4(c)  show  an  example  of  our  relative  NEQ(f)  analysis  of  BP,  MITS  and 
FBP  respectively  with  seven  different  combinations  of  acquisition  parameters. 

x  10*  BP:  nearcw 


Figure  4(a):  Relative  NEQ(f) 
of  BP  for  several  different 
acquisition  combinations 


cycles/mm 

x  MITS  Slabbing5:  nearcw 


Figure  4(b):  Relative  NEQ(f) 
of  MITS  for  several  different 
acquisition  combinations 
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,»  CorFBP  Hamm:  nearcw 


Figure  4(c):  Relative  NEQ(f) 
of  FBP  for  several  different 
acquisition  combinations 


We  found  that  FBP  showed  better  relative  NEQ(f)  results  than  that  of  MITS  for  the  cases 
investigated.  A  reason  is  that  MITS  was  applied  on  traditional  SAA  that  lines  up  the  structures 
on  the  focal  plane  one-dimensionally;  while  FBP  was  applied  on  point-by-point  BP  that  point- 
by-point  corrects  structures  on  the  focal  plane.  FBP  provides  a  trade-off  between  wider  out-of- 
plane  blur  extension  and  reduced  out-of-plane  artifacts.  FBP  with  larger  number  of  projection 
images  is  the  optimal  digital  breast  tomosynthesis  image  reconstruction  method.  Considering  the 
longer  acquisition  time  of  tomosynthesis  sequence  of  49  projection  images,  N=25  number  of 
projection  images  and  TA=50°  without  binning  is  suggested  as  the  optimal  acquisition  mode  for 
FBP  reconstruction. 


The  relative  NEQ(f)  comparison  results  between  MITS  and  FBP  in  this  project  are  only  limited 
to  breast  tomosynthesis  imaging.  The  conclusion  doesn’t  apply  to  other  tomosynthesis 
applications  such  as  chest  imaging  where  MITS  shows  advantages  in  image  quality  and  artifact 
reduction. 


We  have  published  part  of  our  results  in  a  full-length  proceedings  paper  at  SPIE  (See  Appendix 
#4).  We  are  also  working  on  a  Medical  Physics  journal  manuscript  preparation  for  our  relative 
NEQ(f)  analysis  results. 

KEY  RESEARCH  ACCOMPLISHMENTS 

•  Investigated  several  different  reconstruction  algorithms  for  digital  breast  tomosynthesis, 
including  SAA,  NIKL,  BP,  MITS,  FBP  and  MLEM.  Selected  candidate  algorithms  and 
compared  performance  against  each  candidate  algorithm. 

•  Analyzed  impulse  response,  MTF,  and  NPS  by  simulation,  experiments  and 
characterization  to  compare  and  optimize  the  imaging  acquisition  parameters  including 
total  Tomographic- Angle  (TA),  Number  of  projection  images  (N),  and  Reconstruction- 
Slice-Spacing  (RSS). 

•  Investigated  the  importance  of  point-by-point  BP  for  isocentric  motion  in  digital  breast 
tomosynthesis,  especially  for  reconstruction  of  small  structures  such  as 
microcalcifications. 
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•  Investigated  a  new  image  reconstruction  algorithm  of  Gaussian  Frequency  Blending 
(GFB)  to  combine  MITS  and  FBP  together  for  better  breast  tomosynthesis  image 
reconstruction. 

•  Investigated  a  methodology  of  relative  NEQ(f)  analysis  for  comparison  and  optimization 
of  different  imaging  parameters  and  image  reconstruction  algorithms.  Compared  different 
candidate  image  reconstruction  methods  by  relative  NEQ(f)  analysis  and  found  the 
optimal  image  reconstruction  method  of  FBP  with  bigger  number  of  projection  images 
and  wider  tomographic  angular  range. 

•  Ying  Chen  received  a  PhD  in  Biomedical  Engineering  from  Duke  University. 

REPORTABLE  OUTCOMES 

The  following  papers  and  abstract  are  attached  at  appendices  1,  2,  3,  and  4  with  the  same 
numbers.  The  names  of  the  fellow  (Chen)  and  mentor  (Dobbins)  are  boldfaced  for  emphasis. 

1.  Y.  Chen,  J.  Y.  Lo,  and  J.  T.  Dobbins  III.,  “Importance  of  point-by-point  back  projection 
(BP)  correction  for  isocentric  motion  in  digital  breast  tomosynthesis:  Relevance  to 
morphology  of  microcalcifications,”  Med.  Phys.  34(10),  3885-3892  (2007). 

2.  Y.  Chen,  J.  Y.  Lo,  J.  T.  Dobbins  III,  “A  Comparison  between  Traditional  Shift-and  add 
(SAA)  and  Point-by-point  Back  Projection  (BP)  -  Relevance  to  Morphology  of 
Microcalcifications  for  Isocentric  Motion  in  Digit  Breast  Tomosynthesis  (DBT),” 
Proceedings  of  IEEE-Bioinformatics  &  Bioengineering,  2007. 

3.  Y.  Chen,  J.  Y.  Lo,  J.  T.  Dobbins  III,  “Two-dimensional  Shift- And- Add  Algorithm  for 
Digital  Breast  Tomosynthesis  Reconstruction,”  Med.  Phys.  33  (6),  2001  (2006). 

4.  Y.  Chen,  J.  Y.  Lo,  N.  T.  Ranger,  E.  Samei,  J.  T.  Dobbins  III,  “Methodology  of  NEQ(f) 
analysis  for  optimization  and  comparison  of  digital  breast  tomosynthesis  acquisition 
techniques  and  reconstruction  algorithms,”  Proc.  SPIE  6510,  65101-1,  (2007). 

CONCLUSIONS 

This  work  presents  the  comparison  and  optimization  of  several  different  breast  tomosynthesis 
image  reconstruction  algorithms  and  imaging  acquisition  parameters.  We  have  investigated 
several  different  reconstruction  algorithms  including  traditional  SAA,  BP,  NIKL,  MITS,  FBP 
and  MLEM.  We  have  characterized  the  impulse  response,  MTF,  and  NPS  by  simulation, 
experiments,  and  computation  to  compare  and  optimize  candidate  algorithms  and  different 
acquisition  parameters. 

One  important  research  outcome  was  the  comparison  and  optimization  of  algorithms  and 
acquisition  parameters.  We  proposed  a  relative  NEQ(f)  analysis  to  combine  the  signal  and  noise 
performance  together.  The  relative  NEQ(f)  analysis  was  applied  to  DBT  for  the  first  time.  It 
provides  a  means  for  simultaneously  compare  signal  and  noise  properties.  FBP  based  on  point- 
by-point  BP  with  more  projection  numbers  was  founded  to  be  the  best  overall  method  with  the 
highest  relative  NEQ(f)  signal  to  noise  properties. 

Additionally,  we  reported  the  importance  of  point-by-point  BP  for  isocentric  motion  in  digital 
breast  tomosynthesis.  Additionally,  the  importance  of  point-by-point  back  projection  (BP)  for 
isocentric  motion  in  DBT  was  described  and  investigated.  It  improves  the  in-plane  sharpness  of 
structures  such  as  microcalcifications,  which  bear  important  meaning  in  clinical  tasks  for  breast 
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cancer  detection.  A  novel  algorithm  of  Gaussian  Frequency  Blending  (GFB)  to  blend  MITS  and 

FBP  together  was  also  investigated  for  better  image  reconstruction  of  the  breast. 
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Digital  breast  tomosynthesis  is  a  three-dimensional  imaging  technique  that  provides  an  arbitrary  set 
of  reconstruction  planes  in  the  breast  from  a  limited-angle  series  of  projection  images  acquired 
while  the  x-ray  tube  moves.  Traditional  shift-and-add  (SAA)  tomosynthesis  reconstruction  is  a 
common  mathematical  method  to  line  up  each  projection  image  based  on  its  shifting  amount  to 
generate  reconstruction  slices.  With  parallel-path  geometry  of  tube  motion,  the  path  of  the  tube  lies 
in  a  plane  parallel  to  the  plane  of  the  detector.  The  traditional  SAA  algorithm  gives  shift  amounts 
for  each  projection  image  calculated  only  along  the  direction  of  x-ray  tube  movement.  However, 
with  the  partial  isocentric  motion  of  the  x-ray  tube  in  breast  tomosynthesis,  small  objects  such  as 
microcalcifications  appear  blurred  (for  instance,  about  1-4  pixels  in  blur  for  a  microcalcification  in 
a  human  breast)  in  traditional  SAA  images  in  the  direction  perpendicular  to  the  direction  of  tube 
motion.  Some  digital  breast  tomosynthesis  algorithms  reported  in  the  literature  utilize  a  traditional 
one-dimensional  SAA  method  that  is  not  wholly  suitable  for  isocentric  motion.  In  this  paper,  a 
point-by-point  back  projection  (BP)  method  is  described  and  compared  with  traditional  SAA  for  the 
important  clinical  task  of  evaluating  morphology  of  small  objects  such  as  microcalcifications. 
Impulse  responses  at  different  three-dimensional  locations  with  five  different  combinations  of  im¬ 
aging  acquisition  parameters  were  investigated.  Reconstruction  images  of  microcalcifications  in  a 
human  subject  were  also  evaluated.  Results  showed  that  with  traditional  SAA  and  45°  view  angle 
of  tube  movement  with  respect  to  the  detector,  at  the  same  height  above  the  detector,  the  in-plane 
blur  artifacts  were  obvious  for  objects  farther  away  from  x-ray  source.  In  a  human  subject,  the 
appearance  of  calcifications  was  blurred  in  the  direction  orthogonal  to  the  tube  motion  with  tradi¬ 
tional  SAA.  With  point-by-point  BP,  the  appearance  of  calcifications  was  sharper.  The  point-by- 
point  BP  method  demonstrated  improved  rendition  of  microcalcifications  in  the  direction  perpen¬ 
dicular  to  the  tube  motion  direction.  With  wide  angles  or  for  imaging  of  larger  breasts,  this  point- 
by-point  BP  rather  than  the  traditional  SAA  should  also  be  considered  as  the  basis  of  further 
deblurring  algorithms  that  work  in  conjunction  with  the  BP  method.  ©  2007  American  Association 
of  Physicists  in  Medicine.  [DOI:  10.1118/1.2776256] 

Key  words:  mammography,  tomosynthesis,  3D  reconstruction,  shift-and-add  (SAA),  back 
projection,  microcalcifications 


I.  INTRODUCTION 

Breast  cancer  is  the  most  common  cancer  among  women. 
Currently,  mammography  is  the  most  important  and  effica¬ 
cious  tool  for  the  early  detection  of  breast  cancer.1  However, 
limitations  of  mammography  have  been  well  publicized, 
such  as  20%  false  negative  rate,  ’  many  callbacks  from 
screening,  and  low  positive  predictive  value  of  about  15  %- 
34%  from  biopsy.4,5  It  can  be  difficult  for  conventional  two- 
dimensional  (2-D)  mammography  to  distinguish  a  cancer 
from  overlying  breast  tissues. 

Digital  breast  tomosynthesis  is  a  three-dimensional  (3-D) 
imaging  technique  that  provides  an  arbitrary  set  of  recon¬ 


struction  planes  in  the  breast  from  a  limited-angle  series  of 
projection  images  when  the  x-ray  tube  moves.6  There  are  a 
variety  of  tomosynthesis  reconstruction  algorithms,  includ¬ 
ing  the  image-stretching  method  proposed  by  Niklason  and 
colleagues,7  maximum  likelihood  iterative  algorithm  by  Wu 
et  al.f  tuned-aperture  computed  tomography  (TACT)  recon¬ 
struction  methods  developed  by  Webber  and  investigated  by 
Suryanarayanan  et  al.9,10  algebraic  reconstruction 
techniques,11-13  filtered  back  projection  (FBP), 14-20  matrix 
inversion  tomosynthesis  (MITS), 17,21-24  and  Gaussian  fre¬ 
quency  blending  of  MITS  and  FBP.15  Some  of  these  algo¬ 
rithms  depend  on  a  traditional  shift-and-add  (SAA)  method 
that  is  appropriate  for  parallel-path  geometries.  For  example, 
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Niklason  and  colleagues  modified  the  traditional  shift-and- 
add  technique  for  mammography  to  stretch  the  image  along 
the  direction  of  x-ray  tube  motion  to  account  for  the  effects 
of  magnification  variation  with  angle,  but  the  correction  nec¬ 
essary  along  the  direction  perpendicular  to  the  tube  motion 
was  not  taken  into  account.7  Suryanarayanan  et  al.  applied 
Webber’s  TACT  method  to  breast  tomosynthesis 
reconstruction,9,10  and  used  traditional  SAA.  The  MITS  tech¬ 
nique  developed  in  our  laboratory  has  been  investigated  for 

breast  tomosynthesis  using  a  traditional  SAA  algorithm  as 

21 

the  basis  for  subsequent  matrix  inversion  deblurring. 

Traditional  SAA  is  appropriate  for  parallel-path  tube 
movement  when  the  path  of  the  tube  lies  in  a  plane  parallel 
to  the  plane  of  the  detector.6  However,  the  partial  isocentric 
motion  of  the  tube  in  breast  tomosynthesis  causes  a  non¬ 
parallel  motion.  The  path  of  the  tube  does  not  lie  in  a  plane 
that  is  parallel  to  the  detector  plane.  While  the  effects  due  to 
isocentric  motion  are  small  for  most  objects,  the  use  of  SAA 
methods  introduces  morphological  distortions  with  small  ob¬ 
jects  such  as  microcalcifications.  Therefore,  a  simple  SAA 
reconstruction  algorithm  is  not  entirely  suitable  for  breast 
tomosynthesis.  This  issue  of  non-parallel  motion  is  addressed 
in  point-by-point  back  projection  (BP)  methodologies,14  but 
its  impact  has  largely  not  been  evaluated  with  algorithms  that 
rely  simply  on  traditional  SAA  approaches.  This  paper  dem¬ 
onstrates  the  importance  of  point-by-point  corrections  for 
isocentric  motion  in  digital  breast  tomosynthesis  by  examin¬ 
ing  how  the  morphology  of  microcalcification  reconstruc¬ 
tions  changes  relative  to  a  traditional  SAA  method  that  does 
not  employ  point-by-point  corrections. 

In  this  paper,  a  point-by-point  BP  correction  method  is 
described  and  compared  with  traditional  SAA  by  analysis  of 
impulse  response.  Impulse  responses  at  different  3-D  loca¬ 
tions  with  five  different  combinations  of  imaging  acquisition 
parameters  were  investigated.  In  addition,  reconstructed  im¬ 
ages  of  a  calcification  in  a  human  subject  were  evaluated  to 
demonstrate  the  improvement  in  the  morphology  of  micro¬ 
calcifications  associated  with  the  point-by-point  BP  correc¬ 
tion  method. 

II.  METHODS 

II.A.  Breast  tomosynthesis  system 

A  selenium-based  direct  conversion  Siemens  Mammomat 
Novation  DR  prototype  system  was  modified  to  be  used  as 
the  breast  tomosynthesis  acquisition  system.25  The  detector 
area  was  24  cmX30  cm  (2816X3584  pixels),  with  a  pixel 
pitch  of  85  /zm  (different  from  that  used  in  the  clinical  digi¬ 
tal  mammography  system  from  the  same  manufacturer).  Sev¬ 
eral  different  modes  were  provided  to  choose  from  different 
available  projection  numbers,  total  angular  range  and  speeds 
with  exposure  and  readout  cycle  between  0.5  and  0.8  sec¬ 
onds  per  image.  Figure  1  shows  a  diagram  of  the  breast 
tomosynthesis  imaging  system.  During  the  tomosynthesis 
procedure,  the  x-ray  tube  moves  automatically  along  an  arc 
above  the  chest  wall  to  acquire  up  to  49  projection  images.  A 
continuous  x-ray  motion  was  employed.  The  distance  R  from 
the  rotation  center  O  to  the  detector  is  6  cm.  The  range  of  the 


Fig.  1 .  Breast  tomosynthesis  imaging  System.  O  is  the  rotation  center,  R  is 
the  rotation  center  to  detector  distance,  L  is  the  rotation  arm  length,  and  Z  is 
the  height  of  plane  S  above  the  detector. 


actual  view  angles  can  be  up  to  about  45°  with  respect  to  the 
detector  (50°  at  the  rotation  center).  A  compression  paddle  is 
used  to  keep  the  object  still. 

II. B.  Traditional  shift-and-add  algorithm 

The  traditional  SAA  tomosynthesis  reconstruction 
algorithm6,17  is  a  common  mathematical  method  to  line  up 
each  projection  image  based  on  its  relative  shift  to  generate 
reconstruction  slices  at  specified  depths.  When  the  x-ray  tube 
moves,  objects  at  different  heights  above  the  detector  will  be 
projected  onto  the  detector  at  positions  depending  on  the 
relative  heights  of  the  objects. 

In  order  to  reconstruct  slice  images  of  the  breast,  each 
projection  image  should  be  shifted  by  an  amount  appropriate 
for  the  plane  of  reconstruction.  If  the  detector  remains  sta¬ 
tionary  and  the  tube  moves  in  a  plane  that  is  parallel  to  the 
detector  plane,  the  magnification  of  objects  depends  only  on 
the  height  of  the  object.  With  the  traditional  SAA  algorithm 
for  breast  tomosynthesis  reconstruction,  shift  amounts  for 
each  projection  plane  are  calculated  along  the  axis  of  x-ray 
tube  movement.  In  this  paper,  the  shift  amount  was  calcu¬ 
lated  based  on  projected  positions  from  central  points  of  each 
reconstruction  plane.  The  shifted  planes  were  added  together 
to  emphasize  structures  in  the  in-focus  plane  and  blur  out 
structures  in  other  planes. 

As  shown  in  Fig.  1,  plane  S  represents  a  reconstruction 
plane  at  a  height  of  Z  above  the  detector  surface.  When  the 
x-ray  tube  moves,  objects  in  plane  S  will  be  projected  onto 
the  detector  surface.  For  a  specific  projection  image  from 
angle  6 ,  in  order  to  shift  the  projection  image  to  line  up 
structures  from  plane  S ,  the  traditional  SAA  algorithm  uses 
the  shift  amount  calculated  as 

shifty (Z)  =  Pt(Z)  =L  sin  0- -  /  (1) 

L  cos  6  +  (R  -  Z) 

where  L  is  length  of  the  rotation  arm,  and  R  is  the  height  of 
the  rotation  axis  from  the  detector  surface.  One  can  obtain 
the  reconstruction  plane  S  as  the  average  of  all  N  shifted 
projection  images: 
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Fig.  2.  Two-dimensional  arc  path  from  isocentric  tube  motion. 


T(x,y)  =  ^2 
Ni=i 


Ij(x,y)  ®  %y-  shift, (Z)]. 


(2) 


II. C.  Point-by-point  back  projection  algorithm 

With  parallel  tube  motion,  the  traditional  SAA  method 
works  well  as  long  as  the  required  shifting  is  the  same  for  all 
pixels  in  the  reconstruction  plane.  However,  because  of  the 
isocentric  motion  of  the  x-ray  tube,  a  shift  actually  occurs  in 
both  x  and  y  directions  on  each  projection  image,  and  the 
shifting  amount  is  not  uniform  for  all  pixels  on  the  recon¬ 
struction  plane.  Figure  2  shows  the  arc  path  of  motion  when 
the  x-ray  tube  moves  along  the  Y  axis.  Point  A  represents  a 
single  structure  on  a  certain  reconstruction  plane.  Pn  Pp  and 
Pk  are  the  actual  projected  locations  of  point  A  on  the  detec¬ 
tor  with  different  x-ray  tube  locations  of  Tb  Tp  and  Tk.  The 
actual  path  of  projected  locations  follows  a  2-D  arc  rather 
than  a  one-dimensional  line.  Therefore,  to  reconstruct  a 
single  pixel  on  a  reconstruction  plane  at  certain  height  above 
the  detector,  the  shift  amount  should  be  considered  along 
both  x  and  y  directions  for  every  pixel  on  the  reconstruction 
plane. 

With  the  point-by-point  BP  algorithm,  shift  amounts  for 
every  pixel  location  on  each  reconstructed  plane  were  com¬ 
puted,  taking  into  account  the  2-D  arc  projection  location  of 
reconstructed  objects  in  each  plane.  In  Fig.  2,  ( Ax,Ay,Az ) 
represent  coordinates  of  point  A.  ( TxhTyhTzi )  represents  the 
tube  position  along  the  x,  y,  z  axes  when  the  tube  moves  to 
position  Tr  {PxhPyhPzi)  represents  projected  coordinates  of 
points  A  on  the  projection  image.  One  can  calculate  the  two- 
dimensional  shift  amount  as 


P.Xj  =  Txf  - 


(Txi-Ax)(Tzi-  PZj) 
Tzi  -  Az 


Pyt  =  Tyt  - 


C Tyi-Ay)(Tzi-PZi ) 
Tzi  -  Az 


(3) 


Since  Pt  is  located  on  the  detector,  one  can  define  Pzi= 0. 
Thus,  the  above  formula  can  be  simplified  as 


Pxt  =  Txf  + 


C Txi-Ax)Tzi 
Tzi  -  Az 


Pyt  =  Tyl  + 


{Tyi-Ay)Tzi 
Tzi  -  Az 


(4) 

(4) 


The  final  pixel  value  of  point  A  in  the  tomosynthesized  re- 

N 

construction  was  calculated  as  (1  /N)^I(Pt)  where  I(Pj)  is 

i=  1 

the  pixel  value  at  a  given  location  on  the  ith  projection  im¬ 
age,  and  N  is  the  total  number  of  projection  images.  In  this 
paper,  bilinear  interpolation  was  used  to  address  the  issue  of 
partial  pixel  locations.  Computation  times  for  the  point-by- 
point  BP  algorithm  are  roughly  comparable  to  the  SAA 
method.  With  a  computer  of  800  MHz  CPU  and  UNIX  op¬ 
erating  system,  it  takes  less  than  5  min  for  either  the  point- 
by-point  BP  or  traditional  SAA  reconstruction. 


II. D.  Impulse  response  analysis 

A  single  delta  function  was  simulated  by  ray-tracing 
method  as  the  input  impulse  to  investigate  the  sharpness  of 
reconstructed  in-plane  structures  and  to  see  how  the  tradi¬ 
tional  SAA  and  point-by-point  BP  algorithm  differ  from  each 
other.  Two  different  impulse  locations  were  investigated  in 
this  paper:  (1)  an  impulse  that  is  near  the  chest  wall  (20 
pixels  away  from  the  chest  wall)  and  in  a  defined  recon¬ 
structed  plane  (40  mm  above  the  detector  surface);  and  (2) 
an  impulse  that  is  approximately  4  cm  away  from  the  chest 
wall  and  in  a  defined  reconstructed  plane  (40  mm  above  the 
detector  surface).  Parameters  of  the  digital  breast  tomosyn¬ 
thesis  device  described  in  Sec.  II  A  were  used  for  geometries 
of  the  simulation.  Five  different  combinations  of  acquisition 
parameters  including  projection  image  numbers  and  total  an¬ 
gular  range  were  applied:  (1)  13  projections  with  22.7°  view 
angle  range  with  respect  to  the  detector;  (2)  13  projections 
with  45°  view  angle  range  with  respect  to  the  detector;  (3)  25 
projections  with  22.7°  view  angle  range;  (4)  25  projections 
with  45°  view  angle  range  with  respect  to  the  detector;  and 
(5)  49  projections  with  45°  view  angle  range.  For  each  im¬ 
pulse  location  and  combination  of  acquisition  parameters, 
two  datasets  of  projection  images  were  simulated  by  ray¬ 
tracing  method:  (1)  background-only:  only  Hr2  shading  dif¬ 
ference  for  each  pixel  on  projection  images  was  taken  into 
account  (r  is  the  distance  from  the  x-ray  source  to  each  pixel 
location);  and  (2)  impulse- added:  projection  images  with 
simulated  impulse  and  the  Hr 2  shading  difference.  Other 
system  blur  and  noise  issues  are  not  addressed  in  this  paper 
in  order  to  focus  on  the  contribution  of  the  blur  due  to  the 
isocentric  motion.  During  ray-tracing,  if  the  simulated  im¬ 
pulse  was  projected  onto  non-integer  pixel  location  on  the 
detector  surface,  linear  interpolation  of  the  projected  impulse 
among  four  neighboring  pixels  was  performed. 

Traditional  SAA  and  point-by-point  BP  reconstruction  al¬ 
gorithms  were  applied  to  both  impulse-added  and 
background-only  simulated  tomosynthesis  projection  se¬ 
quences.  A  reconstruction  plane  spacing  of  1  mm  was  used. 
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Reconstruction  images  from  background-only  projections 
were  subtracted  from  reconstructions  of  impulse-added  pro¬ 
jections  to  eliminate  background  shading  effects.  Because  of 
the  isocentric  motion,  there  will  be  greater  interpolation  that 
will  affect  the  amplitude  of  the  impulse  response  further 
away  from  the  chest  wall.  In  order  to  compare  impulse  re¬ 
sponses  at  different  locations,  the  impulse  responses  were 
normalized  based  on  the  ideal  condition  when  the  impulse  is 
exactly  located  underneath  the  x-ray  source.  Under  this  ideal 
condition,  the  impulse  was  only  projected  along  the  tube’s 
motion  direction  and  no  interpolation  was  involved  in  the 
direction  orthogonal  to  the  tube’s  motion  direction. 

II.  E.  Human  subject  images 

Human  subject  images  have  been  acquired  on  our  proto¬ 
type  breast  tomosynthesis  system  under  an  IRB -approved 
protocol.  Images  of  one  human  subject  with  notable  calcifi¬ 
cations  were  reconstructed  with  the  traditional  SAA  and 
point-by-point  BP  methods  to  demonstrate  the  effect  of  the 
point-by-point  BP  method  on  reconstructed  calcification 
morphological  appearance.  A  tomosynthesis  sequence  was 
acquired  with  25  projection  images  and  a  total  view  angle  of 
45°  with  respect  to  the  detector  (50°  at  the  rotation  center). 
The  radiographic  technique  for  breast  tomosynthesis  was 
selected  using  technique  optimization  procedures  reported 
previously.26  The  target/filter  for  tomosynthesis  exams 
is  tungsten/rhodium.  The  kilovoltage  setting  is  selected 
to  maximize  a  figure  of  merit  (signal  difference  to  noise  ratio 
squared  per  unit  dose)  for  a  given  breast  thickness  and 
density.  Then  the  mAs  can  be  chosen  to  compare  the  dose 
levels  required  to  achieve  the  same  image  quality  when 
compared  against  the  conventional  molibdium/molibdium  or 
molibdium/rhodium  technique,  or  alternatively  to  compare 
the  image  quality  for  the  same  dose  level.  For  this  specific 
subject,  we  chose  to  do  the  latter.  By  using  the  tomosynthesis 
technique  of  W/Rh  at  28  kVp  (half-value  layer  0.50  mm  Al) 
and  112  mAs  for  this  100%  fatty,  45  mm  breast,  we  main¬ 
tained  the  same  dose  as  the  conventional  left  cranio  caudal 
mammogram. 

III.  RESULTS 

III.A.  Impulse  responses 

Figures  3-5  show  the  impulse  response  results  with  a 
simulated  impulse  located  in  a  defined  reconstruction  plane 
at  40  mm  above  the  detector,  and  approximately  4  cm  away 
from  the  chest  wall.  Figure  3  shows  results  from  simulated 
acquisition  parameters  of  25  projections  and  45°  view  angle 
range  with  respect  to  the  detector.  Figure  4  shows  results 
from  13  projections  and  45°  view  angle  range  with  respect  to 
the  detector.  Figures  5  shows  results  from  25  projections  and 
a  narrower  view  angle  range  of  22.7°  with  respect  to  the 
detector.  On  Figs.  3  through  5,  (a)  and  (c),  give  correspond¬ 
ing  values  in  a  plane  at  the  exact  height  of  the  impulse’s 
location;  (b)  and  (d)  give  the  impulse  response  of  reconstruc¬ 
tion  planes  1  mm  lower  than  the  impulse’s  location;  (a)  and 
(b)  are  results  from  point-by-point  BP;  and  (c)  and  (d)  are 


Fig.  3.  Traditional  SAA  and  point-by-point  BP  impulse  responses:  25  pro¬ 
jection  images  and  45°  view  angle  range  with  respect  to  the  detector,  with 
the  impulse  located  40  mm  above  the  detector  and  about  4  cm  away  from 
the  chest  wall,  (a)  and  (b)  give  the  impulse  response  of  point-by-point  BP; 
(c)  and  (d)  give  the  impulse  response  of  traditional  SAA.  (a)  and  (c)  are 
impulse  responses  in  the  plane  at  the  exact  height  of  the  impulse’s  location; 
(b)  and  (d)  give  corresponding  values  in  a  reconstruction  plane  1  mm  below 
the  impulse’s  location. 


results  from  traditional  SAA.  The  v  and  y  axes  give  the  pixel 
location  on  the  reconstruction  plane,  and  the  plot  displays  the 
normalized  amplitude  of  the  response.  The  X  axis  represents 
the  direction  orthogonal  to  the  tube’s  motion  direction.  The  Y 
axis  represents  the  tube’s  motion  direction.  Only  a  40X40 
pixel  region  close  to  the  impulse  is  shown  for  clarity. 

One  can  see  that  with  traditional  SAA,  when  the  impulse 
is  4  cm  away  from  the  chest  wall  and  the  x-ray  source  moves 


(C)  x  y  (d)  x  y 

Z=40mm  Z=39mm 


Fig.  4.  Traditional  SAA  and  point-by-point  BP  impulse  responses:  13  pro¬ 
jection  images  and  45°  view  angle  range  with  respect  to  the  detector,  with 
the  impulse  located  40  mm  above  the  detector  and  about  4  cm  away  from 
the  chest  wall,  (a)  and  (b)  give  the  impulse  response  of  point-by-point  BP; 
(c)  and  (d)  give  the  impulse  response  of  traditional  SAA.  (a)  and  (c)  are 
impulse  responses  in  the  plane  at  the  exact  height  of  the  impulse’s  location; 
(b)  and  (d)  give  corresponding  values  in  a  reconstruction  plane  1  mm  below 
the  impulse’s  location. 
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Fig.  5.  Traditional  SAA  and  point-by-point  BP  impulse  responses:  25  pro¬ 
jection  images  and  22.7°  view  angle  range  with  respect  to  the  detector,  with 
the  impulse  located  40  mm  above  the  detector  and  about  4  cm  away  from 
the  chest  wall,  (a)  and  (b)  give  the  impulse  response  of  point-by-point  BP; 
(c)  and  (d)  give  the  impulse  response  of  traditional  SAA.  (a)  and  (c)  are 
impulse  responses  in  the  plane  at  the  exact  height  of  the  impulse’s  location; 
(b)  and  (d)  give  corresponding  values  in  a  reconstruction  plane  1  mm  below 
the  impulse’s  location. 

along  a  wider  view  angle  range  of  45°,  the  in-plane  response 
is  noticeably  blurred  and  multiple  peaks  exist  in  a  direction 
that  is  perpendicular  to  the  direction  of  tube  movement  [Figs. 
3(c)  and  4(c)],  reflecting  the  uncorrected  partial  isocentric 
tube  motion.  With  point-by-point  BP,  the  in-plane  response  is 
much  sharper  [Figs.  3(a)  and  4(a)].  When  the  number  of 
projection  images  decreases  to  13,  the  in-plane  response  and 
out-of-plane  blur  become  discrete  [Figs.  4(b)  and  4(d)]  due 
to  limited  projection  numbers.  With  a  narrower  view  angle 
range  of  22.7°,  the  differences  between  traditional  SAA  and 
point-by-point  BP  are  less  obvious  (Fig.  5).  However,  one 
can  still  say  that  the  in-plane  response  of  point-by-point  BP 
is  higher  and  sharper  than  that  of  traditional  SAA. 

Table  I  gives  the  full  width  at  half-maximum  (FWHM) 
measurement  and  full  width  at  a  tenth-maximum  values 
(FWTM)  of  the  in-plane  impulse  responses  along  two  or¬ 
thogonal  directions  when  the  impulse  is  located  at  40  mm 


above  the  detector  surface  and  near  chest  wall.  Table  II  gives 
the  same  measurements  of  FWHM  and  FWTM  when  the 
impulse  is  located  4  cm  away  chest  wall.  When  the  impulse 
is  located  near  the  chest  wall,  there  is  only  a  small  difference 
(<1  pixel  size)  in  FWHM  and  FWTM  values  between  tradi¬ 
tional  SAA  and  point-by-point  BP  for  each  combination  of 
acquisition  parameters.  With  a  narrower  view  angle  range  of 
22.7°,  differences  are  small  too  (less  than  1  pixel  size).  How¬ 
ever,  with  a  wider  view  angle  range  of  45°,  when  the  impulse 
is  located  4  cm  away  the  chest  wall,  major  differences  exist 
along  the  X  axis  (direction  orthogonal  to  tube’s  motion  di¬ 
rection).  Multiple  peaks  and  blurs  appear  along  this  direction 
[Figs.  3(c)  and  4(c)].  Due  to  multiple  peaks  along  the  X  axis, 
FWHM  cannot  adequately  represent  the  real  extension  of 
impulse  response;  for  that  reason,  the  measurement  of 
FWTM  was  also  provided.  The  difference  in  FWTM  values 
is  as  large  as  9  pixels  between  traditional  SAA  and  point-by- 
point  BP.  The  blur  associated  with  tomosynthesis  approaches 
was  estimated  by  the  modulation  transfer  function  (MTF) 
measurement  of  the  middle  projection  (0°)  with  tube 
moving.27  The  estimated  FWHM  and  FWTM  values  from 
blur  including  motion  blur  and  system  blur  were  1.2  and  2.2 
pixels. 


III.B.  Human  subject  images 

Figure  6  shows  regions  of  interest  containing  two  calcifi¬ 
cations  from  images  of  the  left  (L)  breast  of  the  human  sub¬ 
ject.  Reconstructed  structures  by  traditional  SAA  and  point- 
by-point  BP  methods  are  compared.  A  12.75  X  12.75  mm2 
region-of-interest  (ROI)  is  demonstrated.  Compared  with  tra¬ 
ditional  SAA,  edges  of  calcifications  are  sharper  and  better 
focused  with  point-by-point  BP  [Figs.  6(b)  and  6(e)].  A  quan¬ 
titative  measurement  of  the  shape  and  width  of  the  punctate 
calcification  (left-most  calcification)  on  this  ROI  is  given  in 
Table  III.  The  Y  axis  corresponds  to  the  direction  of  x-ray 
tube  motion,  and  the  X  axis  is  the  direction  orthogonal  to  the 
tube  motion  direction.  One  can  see  that,  compared  with  the 
traditional  SAA,  point-by-point  BP  provided  clearer  edge 
shape  and  narrower  width  along  X  axis  that  is  perpendicular 
to  the  x-ray  tube  motion  direction. 


Table  I.  The  impulse  is  located  near  chest  wall.  Full  width  at  half-maximum  and  full  width  at  tenth-maximum  measurements  of  in-plane  impulse  response 
along  two  directions:  tube’s  motion  direction  ( Y  axis)  and  direction  orthogonal  to  tube’s  motion  (A  axis). 


Acquisition  Full  width  at  tenth-maximum  Full  width  at  half-maximum 

parameters  (in  pixels)  (in  pixels) 


Projections 

View 

angle 

range 

(°) 

Traditional 

SAA 

Point-by-point 

BP 

Traditional 

SSA 

Point-by-point 

BP 

X  axis 

Y  axis 

X  axis 

Y  axis 

X  axis 

Y  axis 

X  axis 

Y  axis 

13 

22.7 

2.8 

2.8 

3.4 

2.4 

1.9 

1.2 

1.7 

1.2 

25 

22.7 

2.8 

2.8 

3.5 

2.4 

1.9 

1.2 

1.6 

1.2 

13 

45 

2.7 

2.7 

3.1 

2.1 

1.4 

1.2 

1.3 

1.2 

25 

45 

2.7 

2.8 

3.2 

2.3 

1.4 

1.2 

1.4 

1.2 

49 

45 

2.7 

2.6 

3.2 

2.1 

1.5 

1.2 

1.4 

1.2 
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Table  II.  The  impulse  is  located  4  cm  away  chest  wall.  Full  width  at  half-maximum  and  full  width  at  tenth-maximum  measurements  of  in-plane  impulse 
response  along  two  directions:  tube’s  motion  direction  (7  axis)  and  direction  orthogonal  to  tube’s  motion  (X  axis). 


Acquisition 

Parameters 

Full  Width  at  Tenth  Maximum 
(in  pixels) 

Full  Width  at  Half  Maximum 
(in  pixels) 

View 

angle 

range 

(°) 

Traditional 

SAA 

Point-by-point 

BP 

Traditional  SSA 

Point-by-point 

BP 

Projections 

X  axis  Y  axis  X 

axis  Y  axis 

X  axis  Y  axis 

X 

axis  Y  axis 

13 

22.7 

3.7 

2.8 

2.8 

25 

22.7 

3.6 

2.8 

2.6 

13 

45 

11.8 

2.5 

2.2 

25 

45 

11.7 

2.8 

2.4 

49 

45 

11.4 

2.8 

2.6 

2.5 

2.0 

1.2 

1.2 

1.2 

2.3 

2.0 

1.2 

1.2 

1.2 

2.1 

(~13)a 

Multiple 

1.2 

1.2 

1.2 

peaks 

2.4 

(  13)a 

Multiple 

1.2 

1.2 

1.2 

peaks 

2.1 

(  13)a 

1.2 

1.2 

1.2 

Multiple 

peaks 


aFor  the  45°  scan,  multiple  discrete  peaks  in  the  impulse  response  make  FWHM  not  meaningful  as  a  descriptor  of  impulse  width.  For  these  cases,  the 
approximate  overall  width  of  the  impulse  is  provided. 


Figure  7(a)  is  a  low  dose  middle  projection  image  of  the 
same  human  breast  when  the  x-ray  tube  was  positioned  at  the 
0°  position.  Figures  7(b)  and  7(c)  are  reconstructed  slice  im¬ 
ages  by  traditional  one-dimensional  SAA  and  point-by-point 
BP,  respectively,  at  a  height  of  7.5  mm  above  the  detector. 
The  reconstructed  size  of  the  breast  by  traditional  one¬ 
dimensional  SAA  is  larger  than  that  from  point-by-point  BP 
[Figs.  7(b)  and  7(c)].  This  is  due  to  the  uncorrected  magni¬ 
fication  difference  with  traditional  SAA  for  structures  at  dif¬ 
ferent  locations  in  the  reconstruction  plane.  When  the  struc¬ 
tures  are  located  at  difference  heights  above  the  detector 
surface,  differences  in  magnification  exist.  Structures  at 
higher  locations  above  the  detector  surface  will  be  projected 
onto  the  detector  with  a  larger  size.  Traditional  SAA  does  not 


Fig.  6.  Reconstructed  12.75  X  12.75  mm2  ROI  of  a  human  breast  containing 
calcifications,  Z=18  mm  represents  the  plane  height  closest  to  the  location 
of  the  calcification:  (a)  point-by-point  BP,  Z=16.5  mm;  (b)  point-by-point 
BP,  Z=18  mm;  (c)  point-by-point  BP,  Z=19.5  mm;  (d)  traditional  SAA,  Z 
=  16.5  mm;  (e)  traditional  SAA,  Z=18  mm;  (f)  traditional  SAA,  Z 
=  19.5  mm. 


take  this  magnification  difference  into  account.  After  shift- 
and-add,  structures  at  higher  locations  appear  larger  com¬ 
pared  with  the  same  size  structure  at  lower  location.  The 
reconstructed  structure  size  from  traditional  SAA  reconstruc¬ 
tion  cannot  reflect  the  real  structure  size.  Point-by-point  BP 
correctly  addressed  this  magnification  difference  issue  by 
calculating  shift  amounts  for  every  pixel  location  on  each 
reconstruction  plane.  Therefore,  the  point-by-point  BP  recon¬ 
structed  plane  in  Fig.  7(c)  reflects  the  real  size  and  is  smaller 
than  that  from  traditional  SAA. 

IV.  DISCUSSION 

With  the  traditional  SAA  algorithm  for  digital  breast  to¬ 
mosynthesis  reconstruction,  shift  amounts  for  each  projec¬ 
tion  plane  are  calculated  along  the  direction  of  the  x-ray 
tube’s  movement.  However,  due  to  the  isocentric  motion  of 
x-ray  tube,  shift  amounts  are  not  the  same  for  all  pixels  on  a 
reconstruction  plane  at  certain  height  above  the  detector.  The 
track  of  projected  impulse  locations  actually  occurs  in  two 
dimensions  on  the  detector.  Different  shift  amount  is  re¬ 
quired  by  different  pixel  location  on  a  reconstruction  plane. 
As  a  result,  the  in-plane  structures  are  blurred.  Illustrations 
with  impulse  responses  and  reconstructed  human  subject  im¬ 
ages  demonstrated  that  this  is  an  inherent  problem  of  tradi¬ 
tional  one-dimensional  SAA  with  breast  tomosynthesis.  This 
problem  is  more  obvious  when  the  object  is  farther  away 
from  the  chest  wall  and  higher  above  the  detector.  Thus,  with 

traditional  SAA,  objects  such  as  microcalcifications  appear 

28 

blurred  and  their  apparent  morphology  changes. 

With  point-by-point  BP,  the  artifacts  coming  from  the  iso¬ 
centric  x-ray  tube’s  movement  are  corrected.  The  in-plane 
structures  are  sharper.  While  reconstructions  of  gross 
anatomy  were  adequate  with  either  algorithm,  the  morphol¬ 
ogy  of  small  structures  such  as  microcalcifications  recon¬ 
struction  requires  a  point-by-point  correction.  Morphological 
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Table  III.  Quantitative  measurement  of  the  leftmost  calcification  depicted  in  Fig.  6.  Estimated  accuracy  of 
measurement  is  ±1  pixel. 


Z=  16.5  mm 

Z= 

18  mm 

Z= 19.5  mm 

Shape 

Width 
(in  pixels) 

Shape 

Width 
(in  pixels) 

Shape 

Width 
(in  pixels) 

X 

Y 

X 

Y 

X 

Y 

Traditional 

Curvilinear 

10 

21 

Punctate 

9 

11 

Indistinct 

11 

16 

SAA 

blurred 

Point-by- 

linear 

7 

20 

Punctate 

7 

9 

Curvilinear 

7 

17 

point  BP 

sharp 

artifacts  were  reduced  with  the  point-by-point  correction, 
and  rendition  of  small  objects  such  as  microcalcifications 
were  greatly  improved.  These  results  demonstrate  the  impor¬ 
tance  of  using  a  point-by-point  correction  to  remove  isocen- 
tric  motion  artifacts  in  tomosynthesis  imaging  of  the  breast, 
where  the  morphology  of  microcalcifications  has  an  impor¬ 
tant  bearing  on  clinical  decision  making. 

The  source  of  the  difficulties  with  traditional  SAA  is  the 
variable  magnification  introduced  by  the  partial  isocentric 
motion  of  the  x-ray  tube.  The  magnification  of  objects  in 
different  reconstruction  slices  varies  from  plane-to-plane  as  a 
function  of  the  height  of  slices  above  the  detector.  Even 
within  the  same  reconstruction  plane,  the  magnification  also 
changes  at  different  pixel  locations  due  to  the  partial  isocen¬ 
tric  motion.  Traditional  SAA  does  not  take  this  issue  into 
account.  With  point-by-point  BP,  the  shift  amount  is  calcu¬ 
lated  according  to  the  exact  location  of  the  pixel  in  the  re¬ 
construction  slices,  thereby  addressing  this  issue  of  variable 
magnification. 

None  of  methods  commonly  used  for  breast  tomosynthe¬ 
sis  is  truly  spatially  invariant  due  to  the  partial  isocentric 


(b)  (c) 


Fig.  7.  A  human  breast  demonstrating  a  solitary  calcification:  (a)  low  dose 
middle  (0°)  projection  image  of  the  tomosynthesis  sequence.  The  spectrum 
used  for  the  tomosynthesis  sequence  was  28  kVp  with  W/Rh  target/filter  and 
112.5  mAs  for  a  total  of  25  projection  images  and  45°  view  angle  range 
with  respect  to  the  detector,  (b)  Traditional  SAA  reconstructed  slice  image: 
Z=7.5  mm.  (c)  Point-by-point  BP  reconstructed  slice  image:  Z=7.5  mm. 
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tube  motion,  although  for  structures  close  to  the  chest  wall 
and  near  the  detector,  the  imaging  system  is  approximately 
spatially  invariant.  Therefore,  linear  deblurring  techniques 
will  demonstrate  less  uniform  behavior  in  breast  tomosynthe¬ 
sis  geometries  than  in  parallel-path  tomosynthesis  geom¬ 
etries.  We  found  that  when  the  impulse  is  located  near  the 
chest  wall  or  when  the  total  view  angle  range  is  narrower 
such  as  22.7°  with  respect  to  the  detector,  there  is  no  big 
difference  (<1  pixel  size)  in  FWHM  and  FWTM  values  be¬ 
tween  traditional  SAA  and  point-by-point  BP  for  a  recon¬ 
struction  with  a  height  of  40  mm  above  the  detector.  How¬ 
ever,  with  a  wider  view  angle  range  of  45°,  even  moving  the 
impulse  4  cm  away  chest  wall  shows  a  difference  (about  9 
pixels)  along  the  direction  orthogonal  to  the  direction  of  tube 
motion.  Therefore,  with  a  narrow  angle,  or  for  small  or  thin 
breasts,  the  SAA  algorithm  may  be  tolerated.  However,  with 
a  wide  angle  or  large  breast  size,  the  point-by-point  BP  al¬ 
gorithm  rather  than  traditional  SAA  should  be  used  to  mini¬ 
mize  issues  related  to  isocentric  motion.  Deblurring  algo¬ 
rithms  such  as  FBP  or  MITS,  which  are  an  important 
component  of  high-quality  tomosynthesis  reconstruction, 
should  use  the  point-by-point  BP  method  rather  than  the 
SAA  method  to  generate  the  constituent  images  prior  to  de¬ 
blurring  under  the  same  conditions  of  wide  tube  angle  or 
large  breast  size. 


V.  CONCLUSIONS 

This  work  demonstrates  that  point-by-point  BP  is  an  ef¬ 
fective  method  to  reconstruct  3-D  tomosynthesis  images  of 
the  breast  with  improved  rendition  of  small  structures  such 
as  microcalcifications.  Compared  with  the  traditional  SAA 
algorithm,  the  method  of  point-by-point  BP  takes  into  ac¬ 
count  the  variable  magnification  and  shift  occurring  along 
the  direction  orthogonal  to  tube  movement  due  to  the  isocen¬ 
tric  tube  motion.  Point-by-point  BP  improves  the  sharpness 
and  morphology  of  structures  especially  for  small  objects 
such  as  calcifications.  This  may  prove  helpful  to  radiologists 
in  discriminating  malignant  from  benign  microcalcification 
patterns,  and  thereby  improve  the  accuracy  of  breast  cancer 
detection. 
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Abstract 

Digital  breast  tomosynthesis  (DBT)  is  a  three- 
dimensional  imaging  technique  providing  an  arbitrary 
set  of  reconstruction  planes  in  the  breast  with  limited 
series  of  projection  images.  This  paper  describes  a 
comparison  between  traditional  shift-and-add  (SAA) 
and  point-by-point  back  projection  (BP)  algorithms  by 
impulse  response  and  modulation  transfer  function 
(MTF)  analysis. 

Due  to  the  partial  isocentric  motion  of  the  x-ray 
tube  in  DBT,  structures  such  as  micro  calcifications 
appear  slightly  blurred  in  traditional  shift-and-add 
(SAA)  images  in  the  direction  perpendicular  to  the 
direction  of  tube ’s  motion.  Point-by-point  BP  improved 
rendition  of  micro  calcifications.  The  sharpness  and 
morphology  of  calcifications  were  improved  in  a 
human  subject  images.  A  Filtered  Back  Projection 
(FBP)  deblurring  approach  was  used  to  demonstrate 
deblurred  point-by-point  BP  tomosynthesis  images. 
The  point-by-point  BP  rather  than  traditional  SAA 
should  be  considered  as  the  foundation  of  further 
deblurring  algorithms  for  DBT  reconstruction. 

1.  Introduction 

Breast  cancer  is  second  only  to  lung  cancer  as  the 
leading  cause  of  non-preventable  cancer  death.  Digital 


breast  tomosynthesis  (DBT)  has  been  considered  as  a 
promising  technique  to  improve  early  breast  cancer 
detection.  Compared  with  standard  mammography  two- 
dimensional  imaging  techniques,  digital  breast 
tomosynthesis  (DBT)  three-dimensional  imaging 
methods  improve  conspicuity  of  structures  by  removing 
the  visual  clutter  associated  with  overlying  anatomy 
[1,2,3,4,11].  Several  tomosynthesis  image 
reconstruction  algorithms  have  been  proved  effectively 
in  DBT  image  reconstruction  at  of  yet,  including 
traditional  shift-and-add  (SAA)  [1,2],  Niklason  and 
colleagues’  image  stretching  shift-and-add  [3],  Wu  et 
aT  s  maximum  likelihood  maximization  expectation 
(MLEM)  [4,5],  tuned-aperture  computed  tomography 
(TACT)  [6,7],  filtered  back  projection  [2,16-19], 
matrix  inversion  tomosynthesis  [2,8,12],  algebraic 
reconstruction  techniques  (ART)  [14],  and  Gaussian 
frequency  blending  (GFB)  of  MITS  and  FBP  [9]. 

Traditional  shift-and-add  (SAA)  is  a  common 
mathematical  algorithm  to  calculate  the  shift-amount 
along  x-ray  tube’s  motion  direction  to  line  up  in-plane 
structures  and  suppress  out-of-plane  artifacts  [1,2].  It 
has  been  served  as  the  basis  for  a  few  other  algorithms 
such  as  Niklason  and  colleagues’  image  stretching  SAA 
[3],  TACT  [6,7],  MITS  [2,8],  etc.  With  isocentric 
motion  of  DBT,  the  shift  amount  doesn’t  fall  into  a 
plane  parallel  to  tube’s  motion  plane  [10].  It  actually 
occurs  in  both  x-ray  tube’s  motion  and  the  motion 
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perpendicular  to  tube  motion  direction.  As  a  result, 
structures  such  as  microcalcifications  appear  slightly 
blurred. 

This  paper  compared  the  impulse  responses  and 
modulation  transfer  function  (MTF)  of  traditional  SAA 
and  a  point-by-point  back  projection  (BP)  algorithm. 
Geometries  of  a  Siemens  selenium-based  direct 
conversion  Mammomat  Novation  prototype  system 
were  used  to  simulate  tomosynthesis  projection  images. 
Reconstructed  region  of  interests  (ROI)  of 
microcalcification  of  a  human  subject  were  evaluated. 
A  FBP  deblurring  algorithm  based  on  point-by-point 
BP  was  further  investigated  to  demonstrate  deblurred 
BP  tomosynthesis  images,  as  they  would  normally  be 
portrayed  to  observers. 

Results  showed  that  point-by-point  BP  improved  the 
sharpness  and  morphology  of  microcalcifications 
especially  in  the  direction  perpendicular  to  tube’s 
motion  direction.  The  FBP  deblurring  algorithm  based 
on  point-by-point  BP  provided  good  rendition  of 
detailed  structures,  reduced  out-of-plane  blur,  and 
suppressed  high  frequency  noise.  Point-by-point  BP 
rather  than  traditional  SAA  should  be  considered  as  the 
basis  of  further  deblurring  algorithms  to  improve  the 
sharpness  and  morphology  of  structures  especially  for 
small  objects  such  as  calcifications.  This  may  prove 
helpful  to  radiologists  in  discriminating  malignant  from 
benign  microcalcification  patterns,  and  thereby 
improve  the  accuracy  of  breast  cancer  detection. 

2.  Methods 

Figure  1  shows  a  breast  tomosynthesis  imaging 
acquisition  system  [11]  modified  from  a  selenium- 
based  Siemens  Mammomat  Novation  prototype  system. 
The  x-ray  tube  moves  along  an  arc  above  the  digital 
detector.  The  actually  angular  range  of  the  x-ray  tube 
can  be  up  to  ±25°  at  the  rotation  center  that  is  60  mm 
above  the  detector  surface.  The  pixel  size  is  85  pm.  It 
takes  around  20  seconds  for  a  typical  tomosynthesis 
sequence  acquisition  with  25  projection  images. 

Geometries  of  Siemens  tomosynthesis  imaging 
acquisition  system  were  used  to  simulate  tomosynthesis 
projection  images.  In  this  paper,  a  single  delta  function 
located  at  40  mm  above  the  detector  surface  cover  was 
simulated  as  the  input  impulse  to  evaluate  impulse 
responses  and  modulation  transfer  function  (MTF)  of 
traditional  SAA  and  point-by-point  BP.  Tomosynthesis 
sequence  datasets  of  the  impulse  were  simulated  by 
ray-tracing  method  for  49  projection  images  and  ±25° 


x-ray  tube  angular  range  with  respect  to  the  rotation 
center. 


Figure  1.  Siemens  tomosynthesis  imaging  device 

Two  different  locations  of  simulated  impulses  were 
considered:  1)  the  impulse  was  located  4  cm  away  from 
the  chest  wall  (front  edge  of  the  detector)  and  40  mm 
above  the  detector  surface  cover;  2)  the  impulse  was 
located  near  the  chest  wall  (only  20  pixels  in  distance) 
and  40  mm  above  the  detector  surface  cover.  Shading 
background  effects  of  1/r2  exposure  variation  (r  is  the 
distance  from  the  projection  point  to  the  x-ray  point 
source)  was  taken  into  account.  Projection  images  for 
shading  background  only  and  impulse-added  were 
simulated  separately. 

When  the  x-ray  tube  moves,  objects  at  different 
heights  above  the  detector  will  be  projected  onto  the 
detector  at  different  positions  according  to  the  relative 
heights  of  the  objects  [1].  In  order  to  reconstruct  slice 
images  of  the  breast  at  a  certain  height  above  the 
detector,  each  projection  image  should  be  shifted  by  an 
amount  appropriate  for  the  plane  of  reconstruction. 
Based  on  the  projected  locations  of  the  central  point  on 
each  projection  image,  traditional  shift-and-add  (SAA) 
shifts  each  projection  plane  one-dimensionally  along 
tube’s  motion  direction.  The  shifted  planes  were  added 
together  to  emphasize  structures  in  in-focus  plane  and 
blur  out  structures  in  other  planes  [2].  Point-by-point 
BP  calculates  the  exact  two-dimensional  projected 
location  for  each  pixel  on  reconstructed  planes  at  a 
certain  height  above  the  detector.  All  involved  image 
values  on  contributing  projection  images  were  then 
added  together  and  back  projected  to  generate 
reconstruction  images  [5,  10,  20]. 
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The  simulated  tomosynthesis  datasets  of  projection 
images  were  then  reconstructed  by  traditional  SAA  and 
point-by-point  BP  algorithms.  The  image  size  is  2048 
*2048  with  a  pixel  size  of  85  pm.  A  reconstruction 
slice  thickness  of  5  mm  was  used.  In  order  to  eliminate 
background  shading  effects,  reconstruction  images 
from  background-only  projections  were  subtracted 
from  reconstructions  of  impulse-added  projections.  In¬ 
plane  impulse  responses  on  the  reconstruction  plane 
(40  mm  above  the  detector  surface  cover)  and  out-of- 
plane  responses  (35  mm  above  the  detector  surface 
cover)  were  analyzed  and  compared.  The  impulse 
responses  were  normalized  based  on  the  ideal  condition 
when  the  impulse  is  exactly  located  underneath  the  x- 
ray  source.  Modulation  Transfer  Function  (MTF)  was 
computed  as  normalized  two-dimensional  Fourier 
Transform  of  the  impulse  responses  at  reconstruction 
planes. 


3.  Results 

Figures  2,  3,  4,  and  5  show  impulse  response  and 
MTF  results  of  traditional  SAA  and  point-by-point  BP 
for  simulated  tomosynthesis  projection  images  of  49 
projection  images  and  ±25°  x-ray  tube  angular  range 
(with  respect  to  the  rotation  center).  Figures  2  and  3 
are  results  from  traditional  SAA  and  point-by-point  BP 
respectively  when  the  impulse  was  located  40  mm 
above  the  detector  surface  cover  and  about  4  cm  away 
from  the  chest  wall.  Figures  4  and  5  are  results  when 
the  impulse  was  located  40  mm  above  the  detector 
surface  plate  and  near  the  chest  wall.  In  figures  2 
through  5,  (a)  and  (b)  represent  the  out-of-plane 
responses  in  a  reconstruction  plane  that  is  5  mm  below 
the  impulse’s  location;  (c)  and  (d)  represent  the  in¬ 
plane  responses  in  a  reconstruction  plane  that  is  exactly 
at  the  impulse’s  location.  In  (a)  and  (c),  the  normalized 
amplitude  of  impulse  responses  were  displayed.  The  x 
and  y  axes  represent  the  pixel  location  on  the 
reconstruction  plane,  y  direction  is  the  x-ray  tube’s 
motion  direction.  Only  a  40x40  pixels  region  close  to 
the  impulse  was  shown  for  clarity.  In  (b)  and  (d),  MTF 
curves  were  displayed  with  u  and  v  axes  representing 
spatial  frequencies  conjugate  to  x  and  y. 

When  the  impulse  was  located  4  cm  away  from 
chest  wall,  the  in-plane  response  of  traditional  SAA  is 
noticeably  blurred  in  the  x  direction  that  is 
perpendicular  to  the  direction  of  tube  motion  (figure 
2c).  The  out-of-plane  artifacts  appeared  curved  for 
traditional  SAA,  reflecting  the  uncorrected  partial 
isocentric  tube  motion  (figure  2a).  Point-by-point  BP 


provides  sharper  in-plane  response  without  x  direction 
blur  (figure  3c).  The  out-of-plane  response  changed 
from  curved  to  be  straight  (figure  3  a).  Compared  with 
that  of  traditional  SAA  (figure  2b),  the  MTF 
performance  of  point-by-point  BP  is  much  smoother 
(figure  3b). 

When  the  impulse  was  located  chose  to  chest  wall, 
no  substantial  difference  between  traditional  SAA  and 
point-by-point  BP  was  noticed.  However,  one  can  still 
find  that  the  in-plane  response  of  point-by-point  BP  is  a 
little  sharper  by  giving  a  higher  magnitude  of  the 
impulse  response  (figures  4c  vs  5c). 


Figure  2.  Traditional  SAA  impulse  response  and 
MTF,  impulse  located  about  4  cm  away  from  the 
chest  wall,  (a)  and  (b)  give  the  impulse  response  and 
MTF,  respectively,  at  a  reconstruction  plane  5  mm 
below  the  impulse  location;  (c)  and  (d)  give 
corresponding  values  in  a  plane  at  the  exact  height 
of  the  impulse’s  location. 


Figure  3.  Point-by-point  BP  impulse  response  and 
MTF,  impulse  located  about  4  cm  away  from  the 
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chest  wall,  (a)  and  (b)  give  the  impulse  response  and 
MTF,  respectively,  at  a  reconstruction  plane  5  mm 
below  the  impulse  location;  (c)  and  (d)  give 
corresponding  values  in  a  plane  at  the  exact  height 
of  the  impulse’s  location. 


Figure  4.  Traditional  SAA  impulse  response  and 
MTF,  impulse  located  near  the  chest  wall,  (a)  and 
(b)  give  the  impulse  response  and  MTF, 
respectively,  at  a  reconstruction  plane  5  mm  below 
the  impulse  location;  (c)  and  (d)  give  corresponding 
values  in  a  plane  at  the  exact  height  of  the  impulse’s 
location. 


Figure  5.  Point-by-point  BP  impulse  response  and 
MTF,  impulse  located  near  the  chest  wall,  (a)  and 
(b)  give  the  impulse  response  and  MTF, 
respectively,  at  a  reconstruction  plane  5  mm  below 
the  impulse  location;  (c)  and  (d)  give  corresponding 
values  in  a  plane  at  the  exact  height  of  the  impulse’s 
location. 

The  tomosynthesis  projection  images  of  human 
subjects  have  been  acquired  on  our  prototype  breast 


tomosynthesis  system  (figure  1)  under  an  IRB- 
approved  protocol.  Figure  6  shows  a  reconstructed  7.99 
mm  x  7.99  mm  region  of  interest  (ROI)  of  a  human 
subject  breast  containing  a  round-shaped  calcification. 
The  tomosynthesis  sequence  was  acquired  with  25 
projection  images  and  ±25°  angular  range  at  the 
rotation  center.  The  tomosynthesis  technique  of  W/Rh 
at  28  kVp  and  112  mAs  was  used  to  maintain  the  same 
dose  as  the  conventional  left  cranio  caudal  (LCC) 
mammogram  of  the  human  subject. 


Figure  6.  Reconstructed  7.99  mm  x  7.99  mm  ROI  of 
a  human  breast  containing  a  calcification: 

(a)  point-by-point  BP,  Z=25.5  mm;  (b)  point-by- 
point  BP,  Z=27  mm;  (c)point-by-point  BP,  Z=28.5 
mm; 

(d)  traditional  SAA,  Z=25.5  mm;  (e)  traditional 
SAA,  Z=27  mm;  (f)  traditional  SAA,  Z=28.5  mm. 

In  figure  6,  Z=27  mm  represents  the  in-plane 
reconstruction  of  the  round-shaped  calcification. 
Figures  (a)  and  (d)  represent  the  out-of-plane 
reconstruction  that  is  1.5  mm  below  the  in-plane 
location,  (c)  and  (f)  are  reconstruction  regions  that  is 
1.5  mm  above  the  in-plane  location.  Compared  with 
traditional  SAA  (figure  6e),  edges  of  calcifications  are 
sharper  and  better  focused  with  point-by-point  BP 
(figures  6b).  With  traditional  SAA,  out-of-plane 
structures  were  blurred  and  had  curved  appearance 
(figures  6d  and  6f).  With  2D-SAA,  those  problems 
were  corrected,  and  the  out-of-plane  structures  appear 
as  straight  line  segments  (figures  6a  and  6c). 

Figure  7  shows  a  low  dose  middle  projection  image 
of  the  same  human  breast  when  the  x-ray  tube  was 
positioned  at  the  middle  (0°)  position.  Figures  8  shows 
reconstructed  slice  image  at  a  height  of  16.5  mm  above 
the  detector  surface  cover  by  a  deblurring  FBP  based 
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on  point-by-point  BP.  It  applies  a  deblurring  ramp  filter 
to  projection  images  and  point-by-point  back  projects 
the  filtered  projections  to  generate  reconstruction 
images  [2,9].  A  Hamming  filter  and  a  Gaussian  filter 
was  further  applied  to  control  high  frequency  noise  [9]. 


Figure  7.  Low  dose  middle  (0°)  projection 


Figure  8.  FBP  reconstruction  slice  based  on 
Point-by-point  BP 


Compared  with  low  dose  middle  projection  that  is 
similar  to  a  two-dimensional  standard  mammography 
image,  tomosynthesis  reconstruction  images  provided 
better  appearance  of  the  breast  with  depth  information. 

The  FBP  deblurred  point-by-point  BP  image  (Fig. 
8)  demonstrates  good  rendition  of  breast  anatomy  out 
to  the  skin  line  and  provided  good  rendition  of  detailed 
structures  and  suppressed  high  frequency  noise. 


4.  Discussion 

Due  to  isocentric  motion  of  x-ray  tube  in  DBT,  the 
track  of  projected  locations  of  objects  occurs  in  both 
tube’s  motion  direction  and  a  direction  perpendicular 
to  tube’s  motion  direction.  Traditional  SAA  calculates 
shift  amount  for  each  projection  image  only  along  the 
tube’s  motion  direction.  As  a  result,  the  in-plane 
structures  are  blurred.  Illustrations  with  impulse 
responses,  MTF,  and  reconstructed  human  subject 
images  showed  that  this  is  an  inherent  problem  of 
traditional  SAA  with  breast  tomosynthesis.  Traditional 
SAA  can  tolerate  when  the  impulse  was  close  to  the 
chest  wall.  However,  when  the  object  was  farther  away 
from  the  chest  wall,  problem  existed  with  ±25°  angular 
range.  Therefore,  with  traditional  SAA,  objects  such  as 
microcalcifications  appear  slightly  blurred  in  the 
direction  perpendicular  to  the  direction  of  tube  motion, 
and  their  morphology  changes. 

The  morphology  of  calcifications  coming  from  the 
isocentric  x-ray  tube’s  movement  is  corrected  with 
point-by-point  BP  algorithm.  The  in-plane  structures 
are  sharper.  Isocentric  motion  of  x-ray  tube  in  DBT 
also  introduces  also  introduces  variable  magnification. 
The  magnification  of  objects  at  different  heights  above 
the  detector  varies  from  plane-to-plane.  Even  within 
the  same  reconstruction  plane,  the  magnification  also 
changes  at  different  pixel  locations.  Therefore,  the 
reconstructed  size  of  an  object  by  traditional  SAA 
doesn’t  reflect  the  real  size  of  the  object.  Point-by¬ 
point  calculates  shift  amount  for  each  pixel  according 
to  the  exact  location  of  the  pixel  in  the  reconstruction 
slices,  thereby,  addressing  this  issue  of  variable 
magnification. 

To  evaluate,  whether  deblurring  filters  gave 
reasonable  results  with  point-by-point  BP  (even  though 
not  perfectly  spatially  invariant),  a  FBP  reconstruction 
with  applied  deblurring  filters  to  point-by-point  BP 
images  was  illustrated.  Although  not  completely 
spatially  invariant,  this  approach  subjectively  did  a 
good  job  of  effective  deblurring  with  point-by-point  BP 
and  was  relatively  fast.  Results  with  human  subjects 
demonstrated  subjectively  that  FBP  deblurring  with 
point-by-point  BP  provides  good  rendition  of  detailed 
structures,  reduced  out-of-plane  blur,  and  suppressed 
high  frequency  noise.  Other  deblurring  algorithms, 
such  as  Matrix  Inversion  Tomosynthesis  (MITS)  [2,  8], 
may  also  prove  effective,  though  not  perfectly,  in  this 
spatially  non-invariant  environment,  and  will  be 
investigated  in  the  future. 
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5.  Conclusion 


This  paper  demonstrates  the  importance  of  using  a 
point-by-point  BP  correction  of  isocentric  motion 
artifacts  in  tomosynthesis  imaging  of  the  breast,  where 
the  morphology  of  microcalcifications  bears  a 
significant  meaning  on  clinical  decision  making. 
Compared  with  traditional  SAA,  point-by-point  BP 
improves  the  sharpness  and  rendition  of  morphology  of 
structures  such  as  microcalcifications.  Point-by-point 
BP  provides  better-focused  in-plane  response  and  less 
curvature  artifacts.  Issues  of  variable  magnifications 
associated  with  isocentric  motion  of  DBT  are  also 
addressed  by  point-by-point  BP.  Our  example  of  FBP 
reconstructions  based  on  point-by-point  BP  shows 
better  rendition  of  breast  anatomy  and  good  rendition 
of  detailed  structures  and  suppressed  out-of-plane 
artifacts.  The  point-by-point  BP  rather  than  traditional 
SAA  should  be  considered  as  the  foundation  for  further 
deblurring  algorithms  such  as  FBP. 
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Tomosynthesis  Reconstruction 

Y  Chen.  J  Lo.  J  Baker,  and  J  Dobbins 
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Purpose:  To  investigate  a  two-dimensional  Shift-And-Add  algorithm  for  three-dimensional  digital 
breast  tomosynthsis  reconstruction  to  correct  for  defects  existing  in  the  traditional  Shift-And-Add 
algorithm  that  calculates  only  one-dimensional  shift  amount  along  the  axis  of  x-ray  tube's  motion. 
Method  and  Materials:  With  the  traditional  Shift-And-Add  (SAA)  algorithm  for  breast  tomosynthesis 
reconstruction,  shift  amounts  for  each  projection  plane  are  calculated  only  along  the  axis  of  x-ray  tube's 
movement.  As  a  result,  small  objects  such  as  microcalcifications  appear  slightly  blurred  in  the  direction 
perpendicular  to  the  direction  of  tube  motion.  In  this  project,  a  two-dimensional  SAA  method  was 
developed  to  correct  for  this  phenomenon.  Shift  amounts  for  every  pixel  location  on  each  reconstruction 
plane  were  computed,  taking  into  account  the  2D  arc  projection  location  of  reconstructed  objects  in  each 
plane.  Bilinear  interpolation  was  used  for  partial  pixel  locations.  Impulses  at  different  3-D  locations 
were  simulated  and  a  few  human  subject  tomosynthesis  sequences  were  acquired  for  investigation. 
Results:  Two-dimensional  SAA  demonstrated  the  improvement  in  the  direction  that  is  perpendicular  to 
the  tube  motion  direction.  For  human  subjects,  the  appearance  of  calcifications  from  2D  SAA  was 
sharper  than  traditional  SAA  at  the  direction  orthogonal  to  the  tube  motion  direction.  The  out-of-plane 
artifacts  of  calcifications  changed  from  curved  to  be  straight.  Conclusion:  Two-dimensional  SAA  is  an 
effective  method  to  reconstruct  3D  tomosynthesis  images  of  the  breast.  Compared  with  the  traditional 
SAA,  the  new  method  corrects  for  2D  shift  amounts  coming  from  the  isocentric  tube  motion.  This 
provides  more  accurate  and  reliable  results  compared  with  other  SAA  algorithms.  Conflict  of  Interest: 

Research  sponsored  in  part  by  a  research  grant  from  Siemens  Medical  Solutions.  ©2006  American 
Association  of  Physicists  in  Medicine 
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ABSTRACT 

As  a  new  three-dimensional  imaging  technique,  digital  breast  tomosynthesis  allows  the  reconstruction  of  an  arbitrary 
set  of  planes  in  the  breast  from  a  limited- angle  series  of  projection  images.  Though  several  tomosynthesis  algorithms 
have  been  proposed,  no  complete  optimization  and  comparison  of  different  tomosynthesis  acquisition  techniques  for 
available  methods  has  been  conducted  as  of  yet.  This  paper  represents  a  methodology  of  noise-equivalent  quanta 
NEQ  (f)  analysis  to  optimize  and  compare  the  efficacy  of  tomosynthesis  algorithms  and  imaging  acquisition 
techniques  for  digital  breast  tomosynthesis.  It  combines  the  modulation  transfer  function  (MTF)  of  system  signal 
performance  and  the  noise  power  spectrum  (NPS)  of  noise  characteristics.  It  enables  one  to  evaluate  the 
performance  of  different  acquisition  parameters  and  algorithms  for  comparison  and  optimization  purposes.  An 
example  of  this  methodology  was  evaluated  on  a  selenium-based  direct-conversion  flat-panel  Siemens  Mammomat 
Novation  prototype  system.  An  edge  method  was  used  to  measure  the  presampled  MTF  of  the  detector.  The  MTF 
associated  with  the  reconstruction  algorithm  and  specific  acquisition  technique  was  investigated  by  calculating  the 
Fourier  Transform  of  simulated  impulse  responses.  Flat  field  tomosynthesis  projection  sequences  were  acquired  and 
then  reconstructed.  A  mean-subtracted  NPS  on  the  reconstructed  plane  was  studied  to  remove  fixed  pattern  noise. 
An  example  of  the  application  of  this  methodology  was  illustrated  in  this  paper  using  a  point-by-point  Back 
Projection  correction  (BP)  reconstruction  algorithm  and  an  acquisition  technique  of  25  projections  with  25  degrees 
total  angular  tube  movement. 

Keywords:  mammography,  tomosynthesis,  breast  imaging,  modulation  transfer  function  (MTF),  noise  power 
spectrum  (NPS),  Noise-equivalent  quanta  (NEQ),  impulse  response,  Back  Projection  (BP) 


1.  INTRODUCTION 


Digital  breast  tomosynthesis  is  a  new  technique  to  reconstruct  an  arbitrary  set  of  planes  in  the  breast  from  a  limited- 
angle  series  of  projection  images  when  the  x-ray  source  moves  along  an  arc  above  the  breast.  As  of  yet,  several 
tomosynthesis  algorithms  have  been  proposed  [3’4’10"22],  including  Shift-And-Add  (SAA)  [3,4],  Niklason  and 
colleagues’  publication  in  1997  of  a  tomosynthesis  method  with  the  x-ray  tube  moved  in  an  arc  above  the  stationary 
breast  and  detector  [10],  Wu  et  al.’s  report  in  2003  of  the  maximum  likelihood  iterative  algorithm  (MLEM)  to 
reconstruct  the  three-dimensional  distribution  of  x-ray  attenuation  in  the  breast [11],  and  the  filtered  back  projection 
(FBP)  algorithms  t1213’18-24^  tuned-aperture  computed  tomography  (TACT)  reconstruction  methods  developed  by 
Webber  and  investigated  by  Suryanarayanan  et  afl6\  algebraic  reconstruction  techniques  (ART)  [1‘ 4' 23'24],  filtered 
back  projection  (FBP)  [4’13’18"21],  matrix  inversion  tomosynthesis  (MITS)  [15,26]  and  Gaussian  Frequency  Blending 
(GFB)  algorithm  to  combine  the  MITS  and  FBP  for  better  reconstruction  [22]. 
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The  selection  of  optimal  acquisition  parameters  and  reconstruction  algorithm  plays  an  important  role  in  producing 
better  performance.  However,  there  are  several  factors  involved  in  this  task  and  some  of  them  are  not  individually 
independent.  An  effective  methodology  will  enable  one  to  optimize  selection  of  acquisition  parameters  and  compare 
algorithms  for  various  digital  breast  tomosynthesis  methods.  Recently,  Godfrey  et  al.  proposed  methods  of  MTF  and 
NPS  measurement  to  quantitatively  characterize  the  optimal  acquisition  parameters  selection  for  chest 
tomosynthesis  application  [6,7].  In  2005  and  2006,  we  also  published  the  impulse  response  and  NPS  analysis  for 
different  acquisition  parameters  and  reconstruction  algorithms  [3,4]. 

In  this  paper,  a  methodology  of  NEQ  (f)  analysis  is  proposed  to  optimize  and  compare  the  efficacy  of  tomosynthesis 
algorithms  and  imaging  acquisition  techniques  for  digital  breast  tomosynthesis.  It  combines  the  modulation  transfer 
function  (MTF)  of  system  signal  performance  and  the  noise  power  spectrum  (NPS)  of  noise  characteristics.  It 
enables  one  to  evaluate  the  performance  of  different  acquisition  parameters  and  algorithms  for  comparison  and 
optimization  purposes. 


2.  METHODS 


MTF2(f) 

The  NEQ  (f)  is  defined  as:  NEQ(  f)  = - .  The  modulation  transfer  function  (MTF)  and  the  normalized 

NNPS(f) 

noise  power  spectrum  (NNPS)  are  included  to  evaluate  the  performance  of  the  reconstruction  algorithms  and  the 


acquisition  techniques.  The  normalized  NPS  (NNPS)  is  defined  as:  NNPS(u,  v)  = 


NPS^U'V) 

gain 2 


due  to  the 


logarithmic  transform  in  digital  tomosynthesis.  Hence,  the  NEQ  (f)  can  be  calculated  by: 


NEQ{f)  = 


gairr  -MTF\f) 

npsIJT) 


.  In  digital  breast  tomosynthesis  reconstruction,  the  step  of  logarithmic  transform 


involved,  and  therefore,  the  determination  of  NEQ  (f)  is  more  similar  to  that  of  a  screen-film  system  than  to  a  linear 
digital  detector. 

2.1  MTF  measurement 

The  MTF  measurement  includes  two  parts:  1)  the  system  MTF  of  the  detector;  2)  the  reconstruction  MTF  associated 
with  specific  reconstruction  algorithm  and  acquisition  parameters. 


2.1.1  System  MTF  Measurement 

The  system  MTF  describes  the  measured  MTF  of  the  detector.  An  edge  method  was  applied  at  a  range  of  tube 
angles  to  see  if  there  is  any  difference  with  angle.  In  this  paper,  five  different  angular  positions  of  0°,  ±15°,  and  ±25° 
of  the  X-ray  tube  were  selected  as  the  representative  angles. 

Table  1  shows  the  MW2  technique  we  used  for  the  system  MTF  measurement. 


Technique 

Acquisition 

Mode 

kV 

Spectrum 

Number  of 
views  (N) 

Total  mAs 

mAs  for  a  single 
projection 

MW2 

B0XD11  (slow 
speed) 

28 

W/Rh 

11 

303 

~  27.5 

Table  1.  System  MTF  measurement  technique 
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Figure  1  shows  the  set  up  of  the  system  MTF  measurement  experiment.  A  0.1mm  Pt-Ir  edge  was  placed  in  contact 
with  the  detector  and  oriented  at  a  1°  ~  3°  angle  with  respect  to  the  pixel  array.  Edge  images  were  then  acquired  with 
MW2  technique  and  2mm  Aluminum  filtration.  A  previously  published  routine  [8,9]  was  used  to  analyze  the  edge 
images  in  a  region  around  the  edge  (ROI=256*256  pixels,  pixel  size=85pm)  to  compute  the  presampled  system 
MTF. 


Edge  Device 

i 

1°~3°  rotation 

i _ 

Detector 


Chest  Wall 


Figure  1.  Set  up  of  system  MTF  measurements 

2.1.2  Reconstruction  MTF  measurement 

The  reconstruction  MTF  describes  the  calculated  MTF  associated  with  specific  algorithm  and  acquisition 
parameters.  A  dataset  of  tomosynthesis  projection  images  of  a  simulated  delta  function  with  the  specific  acquisition 
parameters  was  served  as  the  input  signal.  A  ray-tracing  simulation  method  was  used  to  project  the  single  delta 
function  onto  the  detector  to  simulate  the  tomosynthesis  sequence  of  projection  images.  Figure  2  shows  the  ray¬ 
tracing  method  to  simulate  the  projection  images  of  a  delta  function. 


Figure  2.  Ray-tracing  simulation 

In  Figure  2,  a  single  delta  function  at  z  distance  above  the  detector  was  simulated.  The  delta  function  is  projected 
onto  the  detector  at  Pt  when  the  x-ray  tube  moves  to  Tl  location.  Partial  pixel  sharing  was  performed  if  Pt  falls  into 
non-integer  pixel  location.  With  this  ray-tracing  method,  for  a  specific  group  of  acquisition  parameters,  a  dataset  of 
projection  images  from  different  x-ray  tube  locations  can  be  simulated. 
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Projection  images  with  1/r2  background  and  without  background  were  considered,  where  r  is  the  distance  between 
the  x-ray  source  Tt  and  the  pixel  location  Pt  on  each  projection  image.  Datasets  of  simulated  projection  images  with 
1/r2  background  and  without  background  were  separately  reconstructed  by  the  specific  reconstruction  algorithm.  The 
reconstruction  images  without  background  were  then  subtracted  from  the  reconstruction  images  with  background  to 
get  rid  of  effects  associated  with  the  1/r2  background  ray-tracing  simulation. 

The  reconstruction  MTF  can  be  computed  as  the  Fourier  Transform  of  the  point  spread  function  (PSF)  on  the 
reconstruction  plane  where  the  simulated  delta  function  is  located  (z  distance  above  the  detector). 

2.1.3  Total  MTF 

The  total  MTF  involves  two  parts  as  described  in  2.1.1  and  2.1.2.  It  can  be  calculated  as  the  multiplication  of 
measured  system  MTF  (2.1.1)  and  the  reconstruction  MTF  associated  with  reconstruction  algorithm  and  acquisition 
parameters  (2.1.2).  Linear  interpolation  was  performed  to  match  the  frequency  axes  of  the  system  MTF  and 
reconstruction  MTF  for  multiplication. 


2.2  NPS  measurement 

In  order  to  compute  the  noise  power  spectrum  (NPS)  for  specific  acquisition  parameters,  two  phantom  slabs  of 
BR12  (47%  water/  53%  adipose  equivalent)  for  a  total  of  4  cm  thickness  were  put  directly  on  the  detector  cover  to 
mimic  the  breast  tissue  equivalent  attenuation  and  scattered  radiation.  Ten  identical  tomosynthesis  sequences  of  flat 
images  with  the  phantom  slabs  on  the  detector  were  acquired.  The  tomosynthesis  sequences  were  then  reconstructed. 

Mean- subtracted  reconstruction  image  data  were  analyzed  and  compared  on  a  reconstruction  plane  with  the  same 
height  (z  distance  above  the  detector)  as  described  in  the  reconstruction  MTF  measurement  (2.1.2).  The  purpose  of 
studying  mean- subtracted  NPS  is  to  remove  fixed  pattern  noise,  including  structured  noise  and  system  artifacts.  A 
previously  published  method  [1,2]  was  applied  using  64  ROIs  of  size  128x128  to  examine  the  noise  response. 

2.3  Gain  factor  calculation 

In  order  to  normalize  NPS  and  fairly  compare  different  reconstruction  algorithms  and  acquisition  parameters,  the 
gain  factor  should  be  included  in  the  NEQ  (f)  calculation. 

For  an  x-ray  screen-film  imaging  system,  the  gain  factor  can  be  considered  as  gain  =  log10  ( e )  •  y  ,  where  y  is  the 

point  slope  of  the  film  density-log  x-ray  fluence  function  [28,29].  Figure  3  shows  a  typical  screen-film  nonlinear 
response  to  exposure.  Q  represents  the  number  of  incident  quanta  /  mm2.  From  this  nonlinear  response  curve,  one 
can  calculate  the  point  slope  y  for  a  specific  log(Q)  on  this  response  curve. 


With  digital  tomosynthesis,  the  gain  factor  varies  with  different  reconstruction  algorithms  and  acquisition 
techniques.  In  this  paper,  we  focus  on  the  relative  NEQ  (f)  methodology  of  a  single  algorithm  and  will  not  include 
the  gain  factor  calculation.  The  gain  factor  will  be  computed  in  future  work  where  we  will  actually  compare 
different  algorithms  against  one  another. 


Optical  Density 


Log  Relative  Exposure 


Figure  3.  Screen- film  characteristic  curve 
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